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We present results of a simulation study of inelastic hard-disks vibrated in a vertical container. 
An Event-Driven Molecular Dynamics method is developed for studying the onset of convection. 
Varying the relevant parameters (inelasticity, number of layers at rest, intensity of the gravity) we 
are able to obtain a qualitative agreement of our results with recent hydrodynamical predictions. 
Increasing the inelasticity, a first continuous transition from the absence of convection to one connec- 
tive roll is observed, followed by a discontinuous transition to two convective rolls, with hysteretic 
behavior. At fixed inelasticity and increasing gravity, a transition from no convection to one roll can 
be evidenced. If the gravity is further increased, the roll is eventually suppressed. Increasing the 
number of monolayers the system eventually localizes mostly at the bottom of the box: in this case 
multiple convective rolls as well as surface waves appear. We analyze the density and temperature 
fields and study the existence of symmetry breaking in these fields in the direction perpendicular 
to the injection of energy. We also study a binary mixture of grains with different properties (in- 
elasticity or diameters). The effect of changing the properties of one of the components is analyzed, 
together with density, temperature and temperature ratio fields. Finally, the presence of a low- 
fraction of quasi-elastic impurities is shown to determine a sharp transition between convective and 
non- convective steady states. 

PACS: 45.70.Qj 

PACS numbers: 

I. INTRODUCTION 

During the last two decades there has been an upsurge of interest for the physical mechanisms which control the 
behavior of granular media, i.e. systems consisting of a large number of macroscopic grains, such as sand, cereals, 
powders, etc [l| These materials play an important role in many industrial and technological processes and in natural 
phenomena and their handling has developed into a multi-billion-dollar industry. Much experimental and theoretical 
effort has been spent on understanding their behavior under a variety of conditions. Among the most frequently 
studied systems are the so called granular gases |2|, obtained by subjecting to an external driving force an assembly of 
grains, so that their behavior resembles that of a molecular gas. However, important differences between the micro- 
world and the macro-world, i.e. between the atomic scale and the millimeter scale, render the analogy incomplete, 
so that one cannot infer their properties from the knowledge of the molecular level. This irreducibility stems chiefly 
from the presence of non conserving forces. Striking manifestations of the peculiarity of granular gases are the non 
Maxwellian velocity distributions, the shear instability, the cluster formation to mention just a few. 

A typical experiment aimed to probe the behavior of granular gases consists of a vertical container partially filled 
with spherical particles which are accelerated by a vertically vibrating base |3( . The competition between the dissipa- 
tion of kinetic energy, due to inelastic collisions between the grains, and the energy provided by the external driving 
force may lead the system to exhibit a variety of non-equilibrium statistically steady states. The phenomenology 
resulting from such a simple experimental setup is incredibly rich and by no means trivial. In addition, by varying the 
control parameters, such as the number of particles, the driving frequency and amplitude, the container aspect-ratio, 
one can observe the crossover from one dynamical phase to another. The exploration of the resulting phase diagram 
has been conducted by means of laboratory experiments and numerical simulations. Several properties have also been 
obtained by granular hydrodynamics and kinetic theory. However, such calculations are difficult as the walls and the 
moving base of the container together with the gravitational field break the translational invariance of the system. 
As a result, one observes density and temperature gradients which render the theoretical analysis highly non trivial. 
In addition, the description of the boundary layers is out of the range of applicability of hydrodynamics. 

In the present work we shall focus on the thermal convection instability in granular systems 0, 0, 0, 0, @, EJ , which 
consists in the appearance of convection rolls due to the competition between temperature gradients and gravit; 
appears to be different from standard convection, which is induced by boundaries and excluded volume effects 
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Thermal convection in granular media has been first observed in 2D simulations 3 , confirmed by 3D experiments 
and analytically investigated in @, Ej . Some numerical investigations can also be found in [(J Q . The underlying 
mechanism is analogous to Rayleigh-Benard convection in classical fluids 0| , with the remarkable difference that in 
a vibrofluidized granular medium the required temperature gradient sets in spontaneously as a consequence of the 
interplay between the collisional dissipation of energy in the bulk of the granular gas and of the power injected by the 
vibrating base. 

Hereafter, we shall consider the onset of convective rolls in mono and bi-disperse two-dimensional vibrated granular 
gases and its similarity with the corresponding phenomenon in ordinary fluids. In particular we shall investigate, by 
means of event-driven Molecular Dynamics simulations, a series of predictions provided by the linear stability analysis 
of hydrodynamic equations 9] for pure granular materials. We shall address the limits of very low and very high 
intensity of the gravity. In addition, we shall exp lore the effects of pol y dis p ersity on thermal convection, studying a 
binary mixture of granular gases [TJ, [^J Uj, HE [lg, [13, EE HO, IH, u3 w ^ n different mechanical properties 
(in particular different coefficients of restitution). We shall characterize the standard hydrodynamic fields, namely 
velocity, density and temperature, both in the presence and in the absence of convective rolls. 

In section II we define the model employed to simulate the grains in the vibrated container, discuss the relevant 
dimensionless control parameter and recall the hydrodynamic predictions. In section III we investigate the pure 
system and propose a qualitative comparison with the hydrodynamic predictions 0, at high and low intensity of 
the gravity. In section IV we study the binary mixture case, measuring velocity, density, and temperature fields. In 
section V we study the effect of gradually adding a quasi-elastic component to a mono-disperse system composed of 
inelastic disks. We observe a sharp transition between a convective and a non convective regime. Finally in section 
VI we present a brief discussion and the conclusions. 

II. THE MODEL 

Let us consider an assembly of grains constrained to move on a vertical rectangular domain, representing the 
container, of dimensions L x x L Z1 and subjected to a gravitational force acting along the negative z direction. The 
grains are idealized as inelastic hard disks of diameter a and are fluidized by the movement of the base, which oscillates 
with frequency v = co/2tt and amplitude A. The collisions of the particles with the side and top walls conserve their 
kinetic energy since the latter are immobile, smooth and perfectly elastic. Instead, the collisions between grains are 
inelastic, and will be represented by means of non constant coefficients of restitution, o^-, which are functions of the 
pre-collisional relative velocity along the direction joining the centers of a pair, V n: according to the formula: 




for V n < v 



{ rij for V n > v 

where i and j indicate the species of the colliding particles, are constants related to the three types of colliding 
pairs, vo = \f9& and g is the gravitational acceleration. At large relative velocities the functions ol^ assume constant 
values , representing the coefficients of restitution, but tend to the elastic value ol^ (0) = 1 for collisions occurring 
with vanishing relative velocity. We have checked the insensitivity of our results on the choice of the value of vo. The 
use of velocity dependent coefficients of restitution has the merit of avoiding anomalous sequences of collisions among 
"collapsed" particles 26]. All the lengths are measured in terms of a. 

The full dynamics consists of a succession of free streaming trajectories, which in the presence of the gravitational 
field have parabolic shapes, and inter-particle collisions or wall-particle collisions. It is therefore very convenient to 
employ the Event Driven Molecular Dynamics. We used the same simulation code of reference |25j |. 

The domain was divided into cells and in each cell time-averaged quantities where computed in order to measure the 
hydrodynamic fields. In such a coarse-grained description we studied the density, p(x,z), and velocity field, v(x,z). 
The local fluctuation of the velocity defines the granular temperature in the cell C XjZ centered around the point of 
coordinates (x, z)\ 

T(x,,) = |(|v-(v) Cx ,JV, . 

where the symbol (-)c x z is an average performed in the cell C x , z - In the measurement procedure we discarded an 
initial transient and focused on steady state properties. 
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From the analysis of the time averaged velocity field we identified the convection rolls without ambiguity by 
measuring the circulation, <£, of the velocity field : 



R 



C(R) 



The value of <!> is obtained, in a steady state, by calculating the sum of the integrals along N pat hs different circular 
paths C(R) of radius R and of length L(C(R)) around the center of the box. The circulation vanishes in the absence 
of convection or in the presence of an even number of convection rolls, and is finite for an odd- number of rolls. 

To characterize the dynamical states of the system, it is convenient to introduce three relevant dimensionless 
parameters IHEj. 

The first parameter is the Froude number, which quantifies the relative importance of the potential energy due 
to the gravity with respect to the energy input at the base. In many previous numerical studies jj, |6| , the base 
is "thermal": a particle colliding with the base is reinjected with a velocity taken from a Gaussian distribution of 
given temperature T5 ase , and thus the Froude number is given by F r = ™® Lz . In hydrodynamical studies @, E|, the 

temperature at the base, To, is a boundary condition and F r = m ^ Lz . In our model, the base is really vibrating, and 
the energy received on average by a colliding particle is mA 2 uj 2 , so that we can choose as definition for F r : 

^ = |% = —, (2) 

where z max = A 2 co 2 /2g is the maximum height reached by a projectile launched vertically with initial velocity 
= Aw. It is worth to comment that F r is employed in the present analysis instead of T = Aw 2 /g, which represents 
the ratio between the acceleration transmitted by the vibrating wall and the gravitational acceleration. In the case of 
a vibrating base, we note that, in contrast with the less realistic case of a thermal wall 0,0 which has been used in 
previous simulations, the mass of the particles does not appear in the definition of the Froude number; this is a relevant 
point when studying mixtures since different masses do not lead to different Froude numbers. Moreover, the presence 
of a realistic vibrating wall, instead of a stochastic hot wall (which are not equivalent, see e.g. [23|), represents a more 
severe test of the conjectured convective scenario. In particular, one can answer the question whether the movement 
of the base breaks the structures and the correlations and renders the convective rolls unstable. 
The second dimensionless parameter relevant for the problem, the Knudsen number 0, |?j 

K = A ((Ti>)) -i = 2 (3) 

V 71 " VK-N layers 

is related to the mean free path and Ni ayers = Ncr/L is the number of filled layers at rest. Finally, in the case of pure 
(monodisperse) systems, we consider the dissipative parameter 

R = 8qK~ 2 = tt(1 - r)N? ayers (4) 

where q = (1 — r)/2 is a measure of the inelasticity of the system and r is the coefficient of restitution. R depends 
both on the inelasticity and on the collision rate, since R — > if either r — > 1 or Ni ayers — » 0. 

It is useful to recall the hydrodynamic predictions concerning the phase behavior of the system. In 0, Ej it has 
presented the following phase diagram (see figure CJ: at fixed F r and K, convection rolls appear with increasing 
inelasticity, i.e. if R overcomes a critical value R c . Such a value, R Cl is an increasing function of the Knudsen number, 
which in turn decreases with the number of particles present in the system. With respect to the Froude number, 
instead R c is a non-monotonic function of F r . As shown in figure ^ a t low F r (i.e. at low gravity or strong shaking) 
R c first decreases, i.e. convection is easier to obtain as the gravity increases. R c however reaches a minimum and then 
increases as the gravity is further increased. In the next section we show the main results of numerical simulations, 
observing good qualitative agreement with predictions from hydrodynamics; however, since vibrating and thermal 
walls can yield different behavior for the same value of the parameters 6], no quantitative test will be made. 



III. CONVECTION IN MONO-DISPERSE SYSTEMS 

We begin by exploring the phase diagram sketched in figure obtained from the theory @,[9], along a plane where 
K is constant and moving along two different paths: a line with constant F r (path A) and a line with constant R 
(path B), respectively. 
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FIG. 1: Phase diagram in the plane (F r ,R), at fixes K, showing the predictions of hydrodynamic theories H,EJ. Convection 
is expected increasing R, e.g. decreasing the restitution coefficient r, as well as changing F r in an adequate interval. At too 
low or too high values of F r (e.g. very low or very high gravity) the system does not reach convection. The two dashed lines 
indicate the paths followed in the numerical simulations to verify the predictions. Path A corresponds to figure Path B to 
figure |5] 




FIG. 2: Pure system: velocity field averaged between t\ — 250 and £2 = 300, in order to show the onset of convection 
increasing R, with F r = 0.36 and K = 0.071, L x = L z = 113, N = 1800, g = 1, A = 0.5, 00 = 50. From left to right, the three 
values of restitution coefficient used are r = 0.9992, r = 0.96 and r = 0.80, corresponding to R = 0.64, R = 32 and R = 159 
respectively. The length of the arrows is normalized so that the maximum velocity v m ax observed is represented by an arrow 
which has a length approximately equal to the distance between two points of the lattice. This convention is used in all figures. 
In these plots the value of v m ax is equal to 14.3, 14.5 and 13.6 (in cm/sec) for the case R = 0.64, R = 32 and = 159 respectively. 



Along path A the effect of increasing R (i.e. of decreasing coefficient of restitution) is to enhance correlations in the 
system and render steeper and steeper the temperature profile in the vertical direction |9|. This leads to the onset 
of convection and eventually to the increase of the number of convection rolls. In figure [2 we display three pictures 
corresponding to three different values of R. At R = 0.64 one observes the absence of convection, whereas a single 
roll appears at R = 32. Finally two rolls are observed for R = 159. The low values of F r (low gravity) and K (low 
density) allow the particles to explore the whole box, i.e. regions far from the bottom as well as regions near the base. 

Figure [3 displays the evolution of the velocity field circulation during the transient and the stationary regimes in 
the case of R = 159: we see that the system first evolves forming a single convection roll before reaching the stationary 
state characterized by the presence of two rolls. From fig we observe that the formation of a roll occurs sharply. 

In figure EI is shown the circulation (defined in eq. QJ) averaged on times between 300 and 1000 and calculated 
for different restitution coefficients, between 0.80 and 0.99. This picture illustrates the presence of two transitions, 
as the restitution coefficient is reduced: there appears to be a continuous transition from absence of convection to 
one convection roll and a discontinuous transition from one roll to two rolls. This second transition shows hysteresis: 
the solid line is obtained performing simulations with an increasing restitution coefficient, i.e. decreasing R, and 
initializing each run with the last configuration of the previous run. The dashed line is instead obtained decreasing 
the restitution coefficient (increasing R). The transition point changes between the two procedures. As a matter of 
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FIG. 3: Pure system: circulation of velocity field <I> as denned in the text versus time for the rightmost system of figure 
(i.e. F r 0.36 and K 0.071, L x = L z = 113, N 1800, g = 1, A = 0.5, u = 50, r = 0.80, corresponding to R = 159). 
It can be seen that the steady state (identified by <3> c± 0, e.g. two opposite convective rolls) is reached after a long transient 
characterized by non-zero circulation (e.g. an odd number of rolls). In particular the regime between t = 150 and t = 250 
corresponds to one convection roll very similar to the one depicted in the central frame of figure e.g. for a lower value of R. 




FIG. 4: Pure system: time averages (between time 300 and 1000) of circulation $ for different values of restitution coefficient 
r, i.e. different values of R, with N = 1800, L x = L y = 113, K = 0.071, F r = 0.36 (g = 1, A = 0.5, u = 50). The solid line 
is obtained decreasing R and initializing each run with the last configuration of the previous run. The dashed line is obtained 
with the same procedure but increasing R. The continuous transition from the absence of convection to one convective roll, as 
well as the discontinuous transition from one roll to two rolls, can be easily recognized. The second transition shows hysteresis. 



fact there is a window of values of R where the system may evolve to one or two rolls, depending on initial conditions. 
Such a hysteretic behavior has been pointed out in |24| , by means of numerical integration of the equations for granular 
hydrodynamics. 

Along path B (R fixed), instead, the system reveals a non-monotonic behavior |8, 9]. In fact, the convection stems 
from the balance between the vertical negative temperature gradients and gravity. Therefore, lowering the intensity 
of gravity (e.g. inclining the box) reduces such a competition and the value of R c at which thermal convection can 



FIG. 5: Pure system: onset and offset of convection (increasing F r ) at R = 3.5 (r = 0.9955) and K — 0.07 fixed, with 
L x = L z = 113, N = 1800, a = 1, A = 0.5, uj = 50. The three values of Froude number are, from left to right, F r = 1, F r = 21 
and F r — 106, corresponding to g = 3, g = 60 and g — 300 respectively. In correspondence of the lowest value of F r convection 
is absent because R < R c , for F r = 21 convection appears because R > R c and finally for F r = 106 convection disappears 
again because R < R c . In these figures the value of v m ax used to normalize the arrows is 14.3, 4.5 and 2.2 (in cm/sec) for the 
case g — 3, g = 60 and g — 300 respectively. 



be observed increases. On the other side, when F r (i.e. gravity) exceeds a certain threshold, there is a crossover to a 
substantially different regime: the grains tend to remain localized near the vibrating plate. Such a spatial arrangement 
determines an increase in density in the lower part of the box and tends to inhibit convection, because the packing 
reduces the mobility of the grains. For this reason when F r becomes larger, the value of the threshold, R c , increases. 
In figure [HI such a non-monotonic behavior is shown: one can see three averaged velocity fields obtained for three 
different values of F r , while keeping R and K fixed. To conclude, the patterns obtained by Molecular Dynamics are 
consistent with the prediction of |£j based on granular hydrodynamics. 
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FIG. 6: Pure system: convection at high F r number (e.g. strong gravity). We used r = 0.991, L x = 375 (L z is irrelevant in 
this case), N = 18000, a = 1, A = 0.5, uj = 50, g = 28, so that K — 0.023 and R — 65. Top: fixed walls boundary conditions; 
bottom: periodic boundary conditions 
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Next, we have tested the prediction of |9( regarding the convection at very large values of the Froude number, i.e. 
when the particles 28] are confined in a region adjacent to the basal plate. In such a case the typical length-scale 
characterizing the onset of convection suggests the presence of a multi-roll configuration. In figure we display the 
velocity field of a system with large aspect ratio L x /L z = 1.5 and a very high Froude number, F r = 22.5. In this 
regime, almost no collisions with the upper wall occur, and L z is therefore essentially irrelevant. Remarkably, five 
convection rolls are present and a wave-like horizontal profile of the density is observed, revealing a horizontal density 
instability (see also 6]). Note that an odd number of rolls can only be obtained with fixed walls boundary conditions, 
while periodic boundary conditions lead always to an even number of rolls, because neighboring cells must rotate in 
opposite directions. 
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FIG. 7: Pure system: convection at high F r number (strong gravity), with the same parameters of figure El (with fixed wall 
boundary conditions): the density (left) and temperature (right) fields in arbitrary units. 

In order to complete our description of the hydrodynamic properties of the system, we have analyzed the granular 
temperature and density fields shown in figure [3 (in the case of fixed lateral walls). Both fields show clearly the 
horizontal inhomogeneity: the density instability and the corresponding temperature instability can be explained 
in terms of the peculiar clustering mechanism characterizing granular fluids. Inelastic collisions, which occur more 
frequently in dense regions, determine a reduction of the local temperature and consequently a drop of the local 
pressure. The ensuing pressure gradient amplifies the density fluctuation and more particles flow toward that region. 

Note that the three regions, where the density is larger, are located in correspondence of the three "valleys" at the 
free surface visible in figure El The temperature field reflects the symmetries of the density field, and one can see 
that density and temperature are anti-correlated: the denser regions are colder. For instance, the three "islands" of 
hot particles correspond to the three lower densities (i.e. the holes left empty by the depression of the free surface, 
see again figure EJ). However, above a certain height (about 100), the horizontal instability disappears and both 
temperature and density decrease together: this is a region with very low number of particles. 

Direct comparison between figures El and [7| clearly shows the strong correlation between the velocity field and 
the temperature and density fields. In the region separating two convective rolls, two opposite phenomena can be 
observed: (i) if the velocity field is directed upward, i.e. toward the free surface, the density is relatively low and 
one observes a bump on the surface; (ii) if the velocity field is directed downward, the bottom plate confines the 
particles and one observes both a depression in the surface and a high density region close to the bottom. Each 
convection roll therefore displays an asymmetric shape with a high-density region coupled to the downward flow and 
a low-density region coupled to the upward flow. We notice that this phenomenology is slightly different from that 
observed in |6[, where a convective regime and a clustering regimes where studied separately: our simulations indicate 
a strong interplay between density and velocity fields. 

IV. BINARY MIXTURES 

In the present section, we shall consider how the convective properties of a granular fluid can be modified by the 
addition of a second species with different physical properties. First, we shall trigger the convection by tuning the 
parameters of one of the two species, in the case of equal concentrations. We also consider the local density ratio, 
that accounts for segregation. 

While non-equipartition can be determined by differences either of masses, sizes or inelasticities, we have mentioned 
in the introduction that the Froude number does not depend on the masses of the particles when the energy is injected 
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FIG. 8: Binary mixture: two different choices of r<2i determine one or two convective cells; (above) the velocity field, 
(bottom) a snapshot of particles positions (crosses and circles are species 1 and 2 respectively). The parameters used are: 
a = 1, L x = L z = 113, Ni = N 2 = 550, g = 1, A = 0.5, oo 50, m 0.96, so that K 0.12 and F r 0.36. The two values 
used for r 2 2 are 0.90, (left) and 0.60 (right). In the two upper plots the value v m ax used to normalize the arrows is 9.5 and 7.2 
(in cm/sec) for the left and right case respectively. 



through a vibrating wall. We shall therefore focus on mixtures of beads with different inelasticities or sizes. 

We shall begin with a mixture consisting of two kinds of grains with different coefficients of restitution. We denote 
with 7*11 the coefficient of restitution relative to collisions among particles both belonging to species 1, T22 relative to 
collisions between particles of species 2 and 7*12 relative to collisions between unlike species. For the sake of simplicity, 
we shall assume Tyi = (rn +^22)/2. Figure [5| shows how the variation of the inelasticity of only one of the two species 
may change the convective structure: in particular we reduce T22 keeping 7*11 fixed. The species 1 has a restitution 
coefficient r\\ = 0.96 equal to the case of figure [21 which displayed a single convection roll. In a first simulation (left 
frame) the species 2 has a restitution coefficient 7*22 = 0.90, and still a single convection roll is observed. In a second 
simulation, with V22 = 0.60, two convection rolls are found. In the same figure we also display two corresponding 
snapshots, in order to give a hint of the coupling between the rolls and the density: as in the previous section, 
regions of high density are coupled to the downward flow of particles. Moreover the snapshots suggest the absence 
of segregation, as expected for dilute systems. The absence of segregation and the visual inspection of convection 
rolls suggest a first simple hypothesis for the behavior of a binary mixture with different restitution coefficients: an 
equimolar mixture behaves similarly to a pure gas with an effective restitution coefficient that is roughly an average 
between those of the two components. In fact in the first case the "average" restitution coefficient is still near the 
value where a single cell was observed in the pure gas, while in the second case this effective coefficient (around 0.78) 
is similar to the one corresponding to two rolls in the pure gas (rightmost frame of figure EJ. 

As far as the temperature ratio field T1/T2 is concerned, we cannot present a conclusive picture. In fact, the tem- 
perature ratio profiles that can be efficiently measured in situations where the systems are horizontally homogeneous ( 
see for instance refs. [l3l[i3ll23) become extremely noisy due to the presence of transversal structures and convection. 
Our results suggest that the temperature ratio is not appreciably affected by convection rolls or at least the variation 
with respect to the case without rolls cannot be resolved due to the statistical uncertainty. We recall that in granular 
gases the breakdown of equipartition is not necessarily related to segregation and can be predicted on the basis of the 
homogeneous Boltzmann equation for dissipative hard spheres (or disks), both in the absence and in the presence of 
external driving forces. 
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FIG. 9: Binary mixture: a mixture of two different sizes {pil<7\ — 0.8 not far from 1) shows differences in the convective 
pattern, with respect of figure|5] The parameters are the same used for that figure. In this plot the value v m ax used to normalize 
the arrows is 7.3 cm/sec. 

Finally, we have studied how to enhance segregation of a mixture, in a case where the global density was higher and 
the sizes of the particles where different. In particular, we have chosen the parameters corresponding to figure E3 but 
substituted half of the grains with grains having larger diameter (02 /<Ji = 0.8). The velocity field is shown in figure[2l 
The field is averaged over very long times, i.e. few hundreds of oscillations of the plate, and the horizontal asymmetry 
observed appears to be stationary, although, in principle, this pattern is expected to be metastable, in the sense that 
an infinitely long average should produce a horizontally periodic field. Therefore, a slight change in the size of half 
of the particles affects the convective field leaving unaltered the number of rolls, but changing their relative size and 
shape. We have analyzed the density fields of the single species n 1 and ri2 and the density ratio 7 = 712/ni, obtaining 
the results in figure EH Segregation is now evident, appearing as deviations of 7 from 1. In particular it can be noted 
that at the bottom there are more particles of species 1 (the larger particles). This results calls for new experiments 
in this dilute "granular gas" regime: in fact it is opposite to well known experimental findings in the densely packed 
regime 3j. Moreover, it appears that stronger segregation is localized in regions of high density. 
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FIG. 10: Binary mixture: Density fields: left m field and right ratio 7 = 712/m, with the parameters of figure 



V. NON EQUIMOLAR MIXTURE 

Even in a granular mixture the relative concentration, \ = ^2/(^1 + N2), is expected to play an important role. 
Not only it affects the ratio of the partial granular temperatures, but also can induce a transition from a homogeneous 
to a convective velocity regime. This has particular relevance, for instance, in industrial processes that require the 
absence or the presence of convection: in the simulations of this section we show how a small percentage of "doping" 
quasi-elastic particles added to an assembly of more inelastic particles can totally suppress convection. 
By employing as an order parameter to detect the transition the circulation <£, of the velocity field, already defined 
in eq. (Q and varying \ we obtained the curve depicted in figure ITT1 In these simulations, we assumed: m = 0.96, 
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FIG. 11: Averaged circulation in a mixture with L x = L y = 113, g = 1, m = 0.96, r 22 = 0.99 92, N± + Af 2 = 2200, 
corresponding to K = 0.06, F r — 0.36. Doping is measured as percentage of quasi-elastic particles, i.e. Nz/(N\ + N2) 



r 22 = 0.9992, N = Ni + 7V 2 = 2200, if = 0.06, F r = 0.36. While the pure system is in a convective regime, a small 
"doping" (the critical value is estimated to be Xcr ~ 0.12) suppresses convection. This critical value changes if the 
relevant parameters of the system (e.g. K or F r ) are changed. For example we have verified that increasing F r the 
value Xcr increases. 
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FIG. 12: Averaged velocity field (between £1 = 250 and t 2 = 300) in a mixture with L x — L y = 113, g — 1, m = 0.96, 
r 22 = 0.99 92, N = 2200, corresponding to K = 0.06, F r = 0.36. Left frame corresponds to the point A in figure ITT1 fa doping 
value x — 0.11), while right frame refers to point B (x = 0.12) 



In figure [T21 we show two different velocity fields observed below and above the critical value of the doping. It is 
quite striking to observe how the system behavior changes qualitatively in correspondence of a tiny variation of the 
doping parameter. The simple picture of an effective restitution coefficient (which in this case should be a weigthed 
average between the two coefficients, using the relative compositions as weights) could be an interpretation of the 
results, as in the previous section. The transition observed in figure could be explained as a transition driven by 
the effective restitution coefficient when the weights are changed. However the validation of such an hypothesis would 
require systematic measurements, beyond the aim of the present investigation. 

To conclude, we have explored the possibility of an hysteretic behavior: the system has been prepared with a doping 
value x slightly higher than Xcr, but with an initial velocity field corresponding to the convective cell at x < Xcr- We 
have studied the temporal evolution of the circulation $(£): the circulation has initially a finite value, but quickly 
decays to zero, indicating that the convection is rapidly washed out. This behavior excludes therefore the possibility 
of hysteresis or long transients. 



VI. CONCLUSIONS 



In this paper, we have studied a two-dimensional granular gas by means of Event-Driven numerical simulation. We 
employed a realistic model regarding both the description of the grains as inelastic hard disks and their interactions 
with the lateral walls and the horizontally oscillating base. The present results are in good qualitative agreement with 
those predicted by the hydrodynamic theory |8|, |9(- I n particular, we observe the same predicted dynamical regimes 
and confirm the non-trivial shape of the phase diagram of |9( , where the effects of inelasticity, number of monolayers 
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and gravity together influence the appearance or disappearance of convection rolls. Within the convective regime, 
the coupling between the various fields (velocity field, density and granular temperature) have been characterized. 
Besides, we have considered the influence of the bidispersity in composition by studying mixtures with different types 
of grains. This is important in industrial applications, where monodispersity is always difficult to be obtained. We 
found that the addition of a small fraction of quasi-elastic component is able to determine the disappearance of the 
convective rolls. Our results suggest that a granular binary mixture of grains with different inelasticities behaves, from 
the point of view of convection, as a monodisperse system with an effective restitution coefficient which is roughly the 
average (weighted with the composition fractions of the two species) of the coefficients of the components. Further 
investigations are necessary to make this statement more precise. 

VII. ACKNOWLEDGEMENTS 

The numerical code used in the present work was developed by Dr. Ciro Cattuto. We wish to thank him for 
supporting the simulation work. 



[1] H.M. Jaeger and S.R. Nagel, Science 255, 1523 (1992); H.M. Jaeger, S.R. Nagel and R.P. Behringer, Phys. Today April 

1996; H.M. Jaeger, S.R. Nagel and R.P. Behringer, Rev. Mod. Phys. 68, 1259 (1996). 

[2] Granular Gases, volume 564 of Lectures Notes in Physics, T. Poschel and S. Luding editors, Berlin Heidelberg, Springer- 

Verlag (2001). 

[3] J. B. Knight, H. M. Jaeger and S. R. Nagel, Phys. Rev. Lett. 70, 3728 (1993). 

[4] R. Ramirez, D. Risso, and P. Cordero, Phys. Rev. Lett. 85, 1230 (2000). 

[5] R. D. Wildman, J. M. Huntley, and D. J. Parker, Phys. Rev. Lett 86, 3304 (2001). 

[6] P. Sunthar and V. Kumaran, Phys. Rev. E 64, 041303 (2001). 

[7] J. Talbot and P. Viot, Phys. Rev. Lett. 89, 064301 (2002); Physica A 317, 672 (2002). 

[8] X. He, B. Meerson, and G. Doolen, Phys. Rev. E 65, 030301 (R) (2002). 

[9] E. Khain and B. Meerson, Phys. Rev. E 67, 021306 (2003). 

[10] S. Chandrasekhar, Hydrodynamic and Hydromagnetic Stability (Dover, New York, 1981). 

[11] J. B. Knight, E. E. Ehrichs, V. Yu. Kuperman, J. K. Flint, H. M. Jaeger, and S. R. Nagel, Phys. Rev. E 54, 5726 (1996). 

[12] Y. Limon Duparcmeur, These de l'universite de Rennes I (1996). 

[13] K. Feitosa and N. Menon, Phys. Rev. Lett. 88, 198301 (2002). 

[14] W. Losert, D.G.W. Cooper, J. Delour, A. Kudrolli and J.P. Gollub, Chaos 9, 682 (1999). 

[15] R.D. Wildman and D.J. Parker, Phys. Rev. Lett. 88, 064301 (2002). 

[16] V. Garzo and J. Dufty, Phys. Rev. E 60 5706 (1999). 

[17] A. Barrat and E. Trizac, Gran. Matter 4, 57 (2002). 

[18] J. M. Montanero and V. Garzo, Gran. Matter 4, 17 (2002). 

[19] R. Clelland and C. M. Hrenya, Phys. Rev. E 65, 031301 (2002). 

[20] U. Marini Bettolo Marconi and A. Puglisi, Phys. Rev. E 65, 051305 (2002) and Phys. Rev. E 66, 011301 (2002). 

[21] R. Pagnani, U. Marini Bettolo Marconi and A. Puglisi, Phys. Rev. E 66, 051304 (2002). 

[22] A. Barrat and E. Trizac, Phys. Rev. E 66, 051303 (2002). 

[23] R. Soto and M. Malek Mansour, Physica A 327, 88 (2003). 

[24] P. Cordero, R. Ramirez and D. Risso, Physica A 327, 82 (2003). 

[25] D. Paolotti, C. Cattuto, U. Marini Bettolo Marconi and A. Puglisi, Granular Matter 5, 75 (2003). 

[26] N. V. Brilliantov, F. Spahn, J. M. Hertzsch, and T. Poschel, Phys. Rev. E 53, 5382 (1996). 

[27] O.R. Walton, Particulate Two-Phase Flow, M.C. Roco editor, Butterworth-Heinemann, Boston (1993), page 884. 

[28] Hong-qiang Wang, Guo-jun Jin, and Yu-qiang Ma, Phys. Rev. E 68, 031301 (2003). 



